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ABSTRACT 

We use the Aquarius simulations to show that the most massive subhalos in galaxy- 
mass dark matter halos in ACDM are grossly inconsistent with the dynamics of the 
brightest Milky Way dwarf spheroidal galaxies. While the best-fitting hosts of the 
dwarf spheroidals all have 12 < l^nax ^ 25kms~^, ACDM simulations predict at least 
ten subhalos with Vmax > 25kms~^. These subhalos are also among the most massive 
at earlier times, and significantly exceed the UV suppression mass back to 2 ^ 10. 
No ACDM-based model of the satellite population of the Milky Way explains this 
result. The problem lies in the satellites' densities: it is straightforward to match the 
observed Milky Way luminosity function, but doing so requires the dwarf spheroidals 
to have dark matter halos that are a factor of ~ 5 more massive than is observed. 
Independent of the difficulty in explaining the absence of these dense, massive subha- 
los, there is a basic tension between the derived properties of the bright Milky Way 
dwarf spheroidals and ACDM expectations. The inferred infall masses of these galax- 
ies are all approximately equal and are much lower than standard ACDM predictions 
for systems with their luminosities. Consequently, their implied star formation effi- 
ciencies span over two orders of magnitude, from 0.2% to 20% of baryons converted 
into stars, in stark contrast with expectations gleaned from more massive galaxies. 
We explore possible solutions to these problems within the context of ACDM and find 
them to be unconvincing. In particular, we use controlled simulations to demonstrate 
that the small stellar masses of the bright dwarf spheroidals make supernova feedback 
an unlikely explanation for their low inferred densities. 
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1 INTRODUCTION 

It has been over ten years since it was first pointed out 
that the number of satellite galaxies of the Milky Way is 
much smaller than the number of dark matter subhalos 
predicted by cold dark matter (CDM)-based simulations of 
Milky Way-mass systems (Klypin et al. 1999; Moore et al. 
1999). These subhalos formed as independent dark matter 
halos and were subsequently accreted onto the Milky Way's 
halo, so each of them is a potential site of galaxy formation; 
that we only see a small number of Milky Way satellites 
therefore requires explanation. This discrepancy may origi- 
nate from galaxies populating only a subset of the predicted 
subhalos, or from modifications of the CDM model on small 
scales such that many of these dark matter subhalos do not 
actually exist. Understanding which (if either) solution is 
correct will have vastly different implications. 
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In this paper, we focus on the former possibility, which 
associates the bright satellites of the Milky Way with the 
most massive subhalos expected from ACDM simulations. 
These models can be further split into two broad classes: 
those that link a satellite's luminosity to the present day 
mass of its dark matter halo (Stoehr et al. 2002; Hayashi 
et al. 2003; Pefiarrubia et al. 2008), and those that as- 
sume satellites form only in the biggest halos defined at 
some earlier epoch, e.g., reionization (Bullock et al. 2000; 
Kravtsov et al. 2004b; Ricotti & Gnedin 2005; Koposov 
et al. 2009; Okamoto & Frenk 2009). Building on our re- 
cent work (Boylan-Kolchin, Bullock, & Kaplinghat 2011b) 
that has demonstrated a puzzling discrepancy between the 
masses of massive simulated subhalos and the dynamics of 
bright Milky Way satellites, we perform a thorough analysis 
of the consistency of models associating bright satellites with 
massive dark matter subhalos. We show that neither of these 
classes reconciles observations with ACDM simulations, ab- 
sent further modifications: massive subhalos, defined either 
now or at any earlier epoch, are generically too dense to be 
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dynamically consistent with the Milky Way satellites. This 
problem relates to the structure of subhalos, making it dis- 
tinct from the standard missing satellites problem, which is 
usually phrased in terms of abundances. 

In fact, galaxy abundance issues are not limited to the 
Milky Way: it has been known for decades that the dark 
matter halo mass function has a steep low-mass slope while 
the faint end of the galaxy luminosity function is flat (Press 
& Schechter 1974; Schechter 1976), signifying a broad need 
for a non-linear mapping between galaxy luminosity and 
halo mass in cold dark matter models (e.g.. White & Rees 
1978). In this sense, the missing satellites problem can be 
interpreted as the low-luminosity manifestation of the faint 
galaxy issue that has been the focus of significant effort in 
galaxy formation theory throughout the modern era. 

In recent years, many groups have shown that a sim- 
ple statistical association between dark matter halos and 
galaxies can be used to reproduce the clustering of galaxy 
light and mass, and its evolution with redshift, using only 
dark matter structures extracted from A'^-body simulations 
(e.g., Kravtsov et al. 2004a; Conroy et al. 2006; Moster 
et al. 2010; Guo et al. 2010). By enforcing the cumula- 
tive abundance of galaxies to be equal to that of halos - 
n(> Afhaio) ~ n{> M),) - one obtains a direct galaxy- 
halo relation that has been tested to scales as small as 
Mhaio ~ 5 X 10^° Mq (Blanton et al. 2008). Abundance 
matching also apparently works well at reproducing the 
luminosity function of satellite galaxies in the Milky Way 
(Busha et al. 2010; Kravtsov 2010; Lunnan et al. 2012), en- 
gendering optimism that the missing satellites problem can 
be solved naturally in ACDM by a Mt(Mhaio) relation that 
falls steeply towards low halo masses. The problem would 
then be reduced to understanding what causes a strong sup- 
pression of galaxy formation in low-mass halos, and there are 
at least two well-motivated processes that likely contribute 
to this inefficiency. 

The first crucial component is thought to be supernova 
feedback. If even a relatively small amount of the energy 
from supernovae can couple to the interstellar medium of a 
galaxy, it can remove gas from shallow gravitational poten- 
tial wells and thereby regulate star formation in low-mass 
systems (Larson 1974; Dekel & Silk 1986). The second pro- 
cess is reionization of the Universe. After the Universe is 
reionized, the resultant UV background affects the gaseous 
content of dark matter halos both by photo-evaporation 
of existing gas and by suppressing subsequent accretion 
through an increase in the temperature of the intergalactic 
medium (Efstathiou 1992; Thoul & Weinberg 1996; Barkana 
& Loeb 1999; Gnedin 2000; Hoeft et al. 2006; Okamoto et al. 
2008). Combined, these effects reduce the baryon fraction of 
halos after reionization in a time-varying manner, such that 
more and more massive halos are affected as time progresses. 
Photo-heating and photo-evaporation therefore may have an 
important effect on galaxy formation in the host halos of 
the Milky Way's dwarf galaxies (Bullock et al. 2000; Benson 
et al. 2002; Somerville 2002; Kravtsov et al. 2004b; Ricotti 
& Gnedin 2005; Madau et al. 2008b; Bovill & Ricotti 2009; 
Busha et al. 2010; Lunnan et al. 2012). 

Semi-analytic models of galaxy formation have shown 
that supernova feedback and reionization should indeed 
strongly limit the formation of stars in low-mass halos since 
z 6 (Benson et al. 2002; Koposov et al. 2009; Muiioz et al. 



2009; Li et al. 2010; Maccio et al. 2010; Guo et al. 2011; 
Font et al. 2011). While the detailed abundance and chem- 
ical composition of MW satellites produced in such models 
depends on the precise combination of reionization and su- 
pernova feedback (Font et al. 2011), it is nevertheless encour- 
aging that such models can reproduce many properties of 
the MW's satellite population using well-established physics. 
Models and simulations including inherently inefficient star 
formation (due to the difficulties in forming molecular gas in 
low-mass halos) have also reproduced dwarf galaxy scalings 
(Robertson & Kravtsov 2008; Tassis et al. 2008; Gnedin et al. 
2009; Kuhlen et al. 2011). Additional feedback processes 
such as cosmic ray pressure may also help suppress star for- 
mation in low mass halos (Wadepuhl & Springel 2011). 

In light of this substantial progress in modeling and suc- 
cesses in comparing observations to simulations, it is tempt- 
ing to consider the ACDM model a success at reproducing 
the main properties of the satellite population of the Milky 
Way. At least one substantial issue remains, however. As 
we showed in a recent paper (Boylan-Kolchin, Bullock, & 
Kaplinghat 2011b, hereafter BBK), the majority of the most 
massive subhalos - the putative hosts of dwarf spheroidals 
- in the highest resolution A/^-body simulations currently 
available are dynamically inconsistent with all of the dSphs: 
massive Milky Way subhalos predicted by ACDM are too 
dense to host the bright MW dSphs. Considering not only 
the mass, but also the structure, of the dSphs' dark matter 
halos provides an additional constraint beyond pure abun- 
dance matching, and shows that while abundance matching 
can reproduce the luminosity function of the Milky Way 
satellites, it does so at the expense of assigning the bright 
satellites to halos that are too dense to actually host the 
MW dSphs. This surprising result presents a new challenge 
for ACDM models on sub-galactic scales. 

In this paper, we present a thorough comparison of the 
bright {Lv > Lq) MW dSphs and simulated MW sub- 
halos, complementing and expanding our results in BBK. 
A substantial difference between BBK and this work is our 
treatment of subhalo mass profiles: rather than relying on 
mass models for the subhalos inferred from their structural 
properties, we use the full particle data from the numeri- 
cal simulations for our comparison. Section 2 describes the 
observational data and numerical simulations upon which 
our analysis is based. We show direct comparisons of sub- 
halo circular velocity profiles with the observed structure of 
Milky Way dSphs in Section 3, demonstrating that there is 
a significant disagreement between the two. In Section 4, we 
use Bayesian inference to determine the masses of the bright 
dSphs' dark matter hosts and show that the results are in- 
consistent with any model relying on an association of the 
brightest satellites with the most massive dwarfs (defined at 
any epoch). Section 5 explores possible reasons for the mis- 
match between the densities of subhalos and the dynamics 
of dwarfs. 



2 SIMULATIONS AND DATA 
2.1 Simulations 

Our ACDM predictions are based on dark matter halos from 
the Aquarius project (Springel et al. 2008), which consists 



© 2012 RAS, MNRAS 000, 1-17 



The Milky Way's bright satellites m ACDM 3 



Table 1. Properties of the Aquarius simulations. Columns: (1) 
Simulation; (2) virial mass; and (3)-(6) number of subhalos at 
2 = within 300 kpc of the halo's center and having VSnfall > 
(20, 30, 40, 50) km respectively. 



Name 


Mvir [Mq] 


N20 




N40 


N50 


Aq-A 


2.19 X 10^2 


105 


33 


15 


6 


Aq-B 


9.54 X 10" 


60 


16 


7 


1 


Aq-C 


1.99 X 10^2 


81 


28 


12 


4 


Aq-D 


2.19 X 10^2 


111 


31 


15 


10 


Aq-E 


1.39 X 10^2 


85 


25 


11 


3 


Aq-F 


1.32 X 10^2 


99 


29 


12 


5 



of six Milky Way-mass dark matter halos (denoted A-F) se- 
lected at z = from a large cosmological simulation^ and 
resimulated at a variety of mass and force resolutions. The 
specific halos selected for re-simulation were chosen at ran- 
dom after applying a mild isolation cut (additionally, the 
selected halos are predicted by a semi-analytic model to 
host late-type galaxies) and are generally typical of ha- 
los with similar mass (Boylan-Kolchin et al. 2010). Only 
halo A was simulated at the highest resolution (level 1), 
in which the particle mass was mp = 1.7 X 10^ Mq and 
the Plummet equivalent gravitational softening length was 
e = 20.5 pc. All halos were simulated at level 2 resolution, 
with fUp = 6.4 X 10^ - 1.4 x lO" Mq and e = 65.8 pc, result- 
ing in approximately 120 million particles within each halo's 
virial radius. These six constitute our sample of simulated 
dark matter halos. 

The masses of the Aquarius halos are (0.95 — 2.2) x 
10^^ M0, a range that reflects the uncertainty in the true 
value of the MW's virial mass and covers almost all recent 
estimates (see Sec. 5.1 for a more detailed discussion of var- 
ious estimates of the mass of the Milky Way's dark matter 
halo). The exact definition of virial mass Mvir itself varies 
among different authors: it is defined to be the mass of a 
sphere, centered on the halo in question, containing an av- 
erage density A times the critical density of the Universe, 
but different authors adopt different values of A. Through- 
out this paper, we use A = Avir, the value derived from the 
spherical top-hat collapse model (Gunn & Gott 1972; Bryan 
& Norman 1998), which results in Avir ~ 94 at 2; = for 
the cosmology used by the Aquarius simulations (see below) . 
A — 200 and A — 200^1rn{z) are two other common values 
used in the literature. 

For each halo, self-bound substructures were identified 
using the SUBFIND algorithm (Springcl et al. 2001) as de- 
scribed in more detail in Springel et al. (2008). Substructures 
can be characterized by their total bound mass Maub, or by a 
characteristic circular velocity Vmax, defined to be the max- 
imum of the circular velocity Vdrc = y/GM{< r)/r. We will 
typically use Kiax (and r-max, defined as Vcirc(7'max) — Vmax) 
rather than Mgub when discussing subhalos because Vmax 
is less dependent on the specific algorithm used to identify 
subhalos and compute their properties. 

Subhalo catalogs were constructed at each time-step for 
which the full particle information was saved, typically 128 



^ This was presented in Gao et al. (2008) as the 'hMS' simula- 
tion, and was resimulated at ~ lOx higher mass resolution as the 
Millennium-II Simulation (Boylan-Kolchin et al. 2009) 



snapshots per halo. The subhalos were linked across snap- 
shots by merger trees, allowing us to explore the full evo- 
lutionary history for each subhalo in addition to its 2: = 
properties. Motivated by abundance matching models, we 
also compute the epoch Zinfaii, defined to be the redshift at 
which a subhalo's mass is maximized (typically just prior 
to infall onto a larger halo), as well as the circular veloc- 
ity Vinfall = Vmax(Zinfall) and maSS Afinfall = Mvir (^infall) at 

that time. In each simulation, we limit our subhalo sam- 
ple to those within 300 kpc of the center of the host halo 
and having Vmax(z = 0) > lOkms"^. Table 1 summarizes 
some important properties of the Aquarius halos and their 
massive subhalos. 

The Aquarius simulations were performed in the con- 
text of a spatially flat WMAPl cosmological model, with a 
matter density of Q.m =0.25, a baryon density of Q,b = 0.045, 
reduced Hubble parameter h = 0.73, linear power spectrum 
normalization erg = 0.9, and a spectral index of the primor- 
dial power spectrum of = 1; this is the same cosmol- 
ogy used in the Millennium and Millennium-II simulations 
(Springel et al. 2005; Boylan-Kolchin et al. 2009). Analy- 
sis of the WMAP7 data indicates that cts = 0.81 ± 0.03 
and Us = 0.967 ± 0.014 (Komatsu et al. 2011, based on 
their "WMAP seven-year mean" values), both of which are 
somewhat lower than the values used in the Aquarius sim- 
ulations. While reductions in these parameters may change 
the properties of halos predicted by ACDM, BBK showed 
that the Via Lactea II (VL-II) simulation (Diemand et al. 
2008), which was run with erg = 0.74 and Us = 0.951, pre- 
dicts very similar structural properties of massive MW satel- 
lites to those in Aquarius. Note that the erg and values 
of VL-II are actually lower than the WMAP7 values, which 
strengthens the impact of this comparison. A preliminary 
analysis of one MW-mass halo run using WMAP7 parame- 
ters shows that the structure of dark matter subhalos in the 
updated cosmology is very similar to that in the Aquarius 
cosmology (Garrison-Kimmel et al., in preparation). Current 
evidence therefore points to our results being independent 
of ~ 10% changes in cosmological parameters, though it is 
certainly desirable to have a large sample of high resolution 
halos simulated in the WMAP7 cosmology in order to make 
the most precise predictions possible. 

2.2 Observational data 

Our primary data for each dwarf are the measured de- 
projected half-light radius {ri/2) and the dynamical mass 
within this radius (Mi/2)- Strigari et al. (2007a), Walker 
et al. (2009) and Wolf et al. (2010) have shown that the dy- 
namical masses of the MW dSphs within a radius compara- 
ble to their stellar extent are well-constrained by kinematic 
data. In particular, while the radial mass proflle is sensitive 
to factors such as the velocity anisotropy, these uncertain- 
ties are minimized at ri/2 (Wolf et al. 2010), resulting in a 
straightforward and accurate (to approximately 20%) esti- 
mate of Ml/2, or equivalently, K:irc(''i/2): 

=3G-i(aL)ri/2, (1) 

Vcirc(rv2) = ^3(a2J. (2) 

The brackets in the above equations refer to luminosity- 
weighted averages. While Eqn. (1) was derived using the 
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spherical Jeans equation, Thomas et al. (2011) have shown 
that this mass estimator accurately reflects the mass as de- 
rived from axisymmetric orbit superposition models as well. 
This result suggests that Eqns. (1) and (2) are also applica- 
ble in the absence of spherical symmetry, a conclusion that 
is also supported by an analysis of Via Lactea II subhalos 
(Rashkov et al. 2012). 

We focus on the bright MW dSphs - those with Lv > 
10^ 1/0 - for several reasons. Primary among them is that 
these systems have the highest quality kinematic data and 
the largest samples of spectroscopically confirmed member 
stars to resolve the dynamics at ri/2- The census of these 
bright dwarfs is also likely complete to the virial radius of 
the Milky Way (~ 300 kpc), with the possible exception of 
yet-undiscovered systems in the plane of the Galactic disk; 
the same can not be said for fainter systems (Koposov et al. 
2008; ToUerud et al. 2008). Finally, these systems all have 
half-light radii that can be accurately resolved with the high- 
est resolution A'^-body simulations presently available. 

The Milky Way contains 10 known dwarf spheroidals 
satisfying our luminosity cut of Lv > 10^ Lq: the 9 clas- 
sical (pre-SDSS) dSphs plus Canes Venatici I, which has a 
IZ-band luminosity comparable to Draco (though it is sig- 
nificantly more spatially extended). As in BBK, we remove 
the Sagittarius dwarf from our sample, as it is in the pro- 
cess of interacting (strongly) with the Galactic disk and is 
likely not an equilibrium system in the same sense as the 
other dSphs. Our final sample therefore contains 9 dwarf 
spheroidals: Fornax, Leo I, Sculptor, Leo II, Sextans, Ca- 
rina, Ursa Minor, Canes Venatici I, and Draco. All of these 
galaxies are known to be dark matter dominated at r-1/2 
(Mateo 1998): Wolf et al. (2010) find that their dynamical 
mass-to-light ratios at rip range from ~ 10 — 300. 

The Large and Small Magellanic Clouds are dwarf ir- 
regular galaxies that are more than an order of magnitude 
brighter than the dwarf spheroidals. The internal dynamics 
of these galaxies indicate that they are also much more mas- 
sive than the dwarf spheroidals: Vcirc(SMC) = 50— 60kms~^ 
(Stanimirovic et al. 2004; Harris & Zaritsky 2006) and 
Kirc(LMC) = 87±5kms-^ (Olsen et al. 2011). Abun- 
dance matching indicates that galaxies with luminosities 
equal to those of the Magellanic Clouds should have Vinfaii ~ 
80 — lOOkms"^ (BBK); this is strongly supported by the 
analysis of Tollcrud ct al. (2011). A conservative estimate 
of subhalos that could host Magellanic Cloud-like galaxies 
is therefore Vinfaii > 60kms~^ and Vmax > 40kms~^. As in 
BBK, subhalos obeying these two criteria will be considered 
Magellanic Cloud analogs for the rest of this work. 



3 COMPARING ACDM SUBHALOS TO 
MILKY WAY SATELLITES 

3.1 A preliminary comparison 

Density and circular velocity profiles of isolated dark mat- 
ter halos are well-described (on average) by Navarro et al. 
(1997, hereafter, NFW) profiles, which are specified by two 
parameters - i.e., virial mass and concentration, or Vmax 
and Tmax- Average dark matter subhalos are also well-fitted 
by NFW profiles inside of their tidal radii, though recent 
work has shown that the 3-parameter Einasto (1965) profile 
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Figure 1. Observed Vcirc values of the nine bright dSphs 
(symbols, with sizes proportional to log Ly), along with ro- 
tation curves corresponding to NFW subhalos with Vmax = 
(12, 18, 24, 40)kms~^. The shading indicates the la scatter in 
'fmax at fixed Vmax taken from the Aquarius simulations. All of 
the bright dSphs are consistent with subhalos having Vmax < 
24kms~^, and most require Vmax J; 18kms~^. Only Draco, the 
least luminous dSph in our sample, is consistent (within 2cr) with 
a massive CDM subhalo of 40kms-l at z = 0. 



provides a somewhat better match to the profiles of both 
simulated halos (Navarro et al. 2004; Merritt et al. 2006; 
Gao et al. 2008; Ludlow et al. 2011) and subhalos (Springel 
et al. 2008) even when fixing the Einasto shape parameter 
(thereby comparing models with two free parameters each). 
To connect this work to the analysis of BBK, Figure 1 com- 
pares the measured values of Vcirc(7'i/2) for the nine bright 
MW dSphs to a set of dark matter subhalo rotation curves 
based on NFW fits to the Aquarius subhalos; the shaded 
bands show the 1 a scatter from the simulations in r-max at 
fixed Vmax. More detailed modeling of subhalos' density pro- 
files will be presented in subsequent sections. 

It is immediately apparent that all of the bright dSphs 
are consistent with NFW subhalos of Vmax ~ 12 — 24 km s~^, 
and only one dwarf (Draco) is consistent with Knax > 
24kms^^. Note that the size of the data points is pro- 
portional to galaxy luminosity, and no obvious trend exists 
between L and Vcirc(''i/2) or Vmax (see also Strigari et al. 
2008). Two of the three least luminous dwarfs, Draco and 
Ursa Minor, are consistent with the most massive hosts, 
while the three most luminous dwarfs (Fornax, Leo I, and 
Sculptor) are consistent with hosts of intermediate mass 
CVnax ~ 18 — 20kms~^). Each of the Aquarius simulations 
contains between 10 and 24 subhalos with Vmax > 25kms~^, 
almost all of which are insufficiently massive to qualify as 
Magellanic Cloud analogs, indicating that models populat- 
ing the most massive redshift zero subhalos with the bright- 
est MW dwarfs will fail. 
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Figure 2. Left panel: circular velocity profiles at redshift zero for subiialos of the Aquarius B halo (top; Afvir = 9-5 X 10^^ ^^'^ E 

halo (bottom; M^ir = 1-4 X 10^^ Mq) that have Vinfall > 30kms^^ and Vmn^iz = 0) > lOkms^^ (excluding MC candidates). Measured 
Vcirc(*'i/2) values for the MW dSphs are plotted as data points with error bars. Each subsequent panel shows redshift zero rotation 
curves for subhalos from the left panel with the ten highest values of Vmax(z = 0) {second panel), Vinfall {third panel), or Knax(z = 10) 
{fourth panel). In none of the three scenarios are the most massive subhalos dynamically consistent with the bright MW dSphs: there 
are always several subhalos more massive than all of the MW dSphs. (Analogous results are found for the other four halos.) 



3.2 Assessing the consistency of massive ACDM 
subhalos with bright Milky Way satellites 

The analysis in Sec. 3.1, based on the assumption that sub- 
halos obey NFW profiles, is similar to the analysis presented 
in BBK. On a case-by-case basis, however, it is possible that 
subhalos may deviate noticeably from NFW profiles. Conse- 
quently, the remainder of our analysis is based on properties 
of subhalos computed directly from the raw particle data. We 
employ a correction that takes into account the unphysical 
modification of the density structure of simulated subha- 
los due to force softening; this procedure is detailed in Ap- 
pendix A. We note, however, that our results do not change 
qualitatively if we neglect the softening correction (see Ap- 
pendix A and Table Al). By using the particle data directly, 
we remove any uncertainties originating from assumptions 
about the shape of the subhalos' density profiles. 

The consistency between massive ACDM subhalos and 
the bright dSphs of the MW is assessed in Figure 2. As there 
is strong theoretical motivation to believe it is Vinfaii rather 
than Vmax(z = 0) that correlates with galaxy luminosity, 
we focus on the most massive subhalos in terms of V^nfaii - 
those with Vinfaii > 30kms~^. We remove from this group 
all subhalos that are Magellanic Cloud analogs according 
to the criteria given at the end of Sec. 2.2. The left-hand 
panels of the figure show circular velocity profiles of the 



remaining massive subhalos in two of the Aquarius halos, 
Aq-B (upper panels; Mvir = 9.5 x 10^^ Mq, the lowest of 
the Aquarius suite) and Aq-E (lower panels; Mvir ~ 1.39 x 
10^^ Mq). Subsequent panels show the ten most massive of 
these subhalos as measured at z — (second column), z = 
Zinfaii (third column), and z = 10 (forth column). 

The most massive subhalos in terms of Vintaii span a 
range of profiles a,t z — 0, as the left panel of Fig. 2 shows. 
For each halo, some of these massive subhalos are consistent 
with the observed data while others are not. Focusing on 
the most massive subhalos at the present day (second panels 
from left), we see that these halos are markedly inconsistent 
with the dSphs, re-enforcing the results of Sec. 3.1. However, 
most subhalos that are massive at z = were also massive in 
the past, a point that is emphasized in the two right panels 
of the figure: the bright MW dSphs are also inconsistent with 
either the most massive subhalos in terms of Mnfaii or those 
defined by their mass at z = 10 (a possible proxy for the 
mass at reionization). Even for Aq-B, the lowest mass host 
halo in the sample, four of the ten most massive subhalos 
are more massive than any of the dSphs, independent of the 
definition of subhalo mass. 

The agreement between MW dSphs and massive subha- 
los is even worse for the other five Aquarius halos. In Fig. 3, 
we compare the redshift zero rotation curves of subhalos 
from each of the six Aquarius halos to the observed values 



© 2012 RAS, MNRAS 000, 1-17 



6 M. Boylan-Kolchin, J. S. Bullock and M. Kaplinghat 




0.1 0.3 0.6 1.0 0.1 0.3 0.6 1.0 0.1 0.3 0.6 1.0 

r [kpc] r [kpc] r [kpc] 



Figure 3. Rotation curves for all subhalos with Vinfall > 30kms^^ and Vmax > lOkms^^, after excluding Magellanic Cloud analogs, in 
each of the six Aquarius simulations (top row, from left to right: A, B, C; bottom row: D, E, F). Subhalos that are at least 2(t denser 
than every bright MW dwarf spheroidal are plotted with solid curves, while the remaining subhalos are plotted as dotted curves. Data 
points with errors show measured Vcirc values for the bright MW dSphs. Not only does each halo have several subhalos that are too 
dense to host any of the dSphs, each halo also has several massive subhalos (nominally capable of forming stars) with Vcirc comparable 
to the MW dSphs that have no bright counterpart in the MW. In total, between 7 and 22 of these massive subhalos are unaccounted for 
in each halo. 



of VciTc{ri/2) for the bright Milky Way dwarf spheroidals. 
As in Fig. 2, we plot only halos with 'Knfaii > 30kms~^ 
and Knax(z = 0) > lOkms"^. Subhalos that are at least 2a 
more massive than every dwarf (at ri/2) are plotted as solid 
curves; these are the "massive failures" discussed in BBK, 
and each halo has at least four such subhalos. Fig. 3 shows 
that each halo has several other subhalos with l^nfaii > 30 
that are unaccounted for as well: for example, halo B has 
three subhalos that are not massive failures by our defini- 
tion but that are inconsistent at 2 cr with every dwarf except 
Draco. Even ignoring the subhalos that are completely un- 
accounted for (and are yet more massive than all of the MW 
dSphs), the remaining massive subhalos do not resemble the 
bright MW dSph population. 



3.3 High redshift progenitors of massive subhalos 

To investigate the possible impact of reionization on our re- 
sults, we show the evolution of the progenitors of all subhalos 
with Vinfall > 30kms~^ in Figure 4. The solid curve show 
the median M{z), while the shaded region contains 68% of 
the distribution, centered on the median, at each redshift. 



For comparison, we also show Tvir{z) = 10* K (the tempera- 
ture at which primordial gas can cool via atomic transitions) 
and 10^ K (dashed lines), as well as the mass Mc{z) below 
which at least half of a halo's baryons have been removed 
by photo-heating from the UV background (Okamoto et al. 
2008). Subhalos with T/infaii > 30kms~^ lie above Mc and 
Tvir ~ 10* K at all redshifts plotted, indicating that they are 
too massive for photo-ionization feedback to significantly al- 
ter their gas content and thereby inhibit galaxy formation. 

Figure 5 focuses on the z = 6 properties of these sub- 
halos. It shows the distribution of halo masses a.t z = 6 
for "massive failures" (open histogram) and the remaining 
subhalos (filled histogram), which are possible hosts of the 
MW dSphs. The massive failures are more massive at 2 = 6, 
on average, than the potentially luminous subhalos. This 
further emphasizes that reionization is not a plausible ex- 
planation of why the massive failures do not have stars: the 
typical massive failure is a factor of ten more massive than 
the UV suppression threshold at a = 6. Implications of this 
result will be discussed in Boylan-Kolchin et al. (in prepa- 
ration) . 

In a series of recent papers, Broderick, Chang, and 
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Figure 4. The median mass of 2 = subhalos having Vinfall > 
SOkms^"*^ (excluding MC analogs) as function of redshift (solid 
curve), along with the 68% confidence range, symmetric about 
the median (shaded region). The hatched region marked "UV- 
suppressed" shows where halos are expected to have lost at 
least 50% of their baryonic mass owing to the UV background 
(Okamoto et al. 2008). Subhalos at z = having Vinfaii > 
30kms~^ are more massive than the photo-suppression scale at 
all redshifts. 



Pfrommer have postulated that the thermal history of the 
IGM at late times {z <2 — 3) could differ substantially from 
standard reionization models owing to a large contribution 
from TeV blazars. This modification relies on plasma in- 
stabilities dissipating energy from TeV blazars in the IGM, 
heating it to a temperature that is a factor of ~ 3 — 10 higher 
than in the case of pure photo-ionization heating. Such heat- 
ing would effectively increase the value of Mc{z) for 2 < 2—3, 
suppressing the stellar content of more massive halos. How- 
ever, as Figure 4 shows, all halos with VinCsii > 30kms~^ 
should have been able to form stars before this epoch, i.e., 
TeV blazar heating happens too late to kill off star formation 
in the Milky Way's massive subhalos (recall that we find 16- 
33 subhalos with l^nfaii > 30kms"^ per halo). While TeV 
blazar heating therefore may help reduce the counts of void 
galaxies (which form later than MW subhalos) and suppress 
the star formation at low redshifts in progenitors of MW sub- 
halos, it does not seem capable of explaining the structure 
and abundance of massive MW satellites. 



4 DARK MATTER MASSES OF THE MILKY 
WAY DWARFS 

The results of Sec. 3 show that the brightest Milky Way 
dwarfs do not inhabit the most massive dark matter subha- 
los from numerical simulations. We can also use the simula- 
tions to compute the properties of subhalos that are consis- 
tent with the dynamics of the bright dSphs. These calcula- 




Figure 5. The distribution of masses at z = 6 for z = subhalos 
with Vinfaii > 30kms~^ (excluding MC candidates). The open 
black histogram shows the "massive failures" (subhalos that are 
2 a more dense than all of the M W's bright dSphs) , while the filled 
gray histogram shows the remaining subhalos, which we deem to 
be potentially luminous satellites. Even at 2 = 6, the massive 
failures are typically more massive than subhalos consistent with 
the bright dSphs. They are also all more massive than character- 
istic scale below which the UV background significantly affects 
the baryon content of halos (hatched region). At z = 6, this char- 
acteristic mass is comparable to Tvij = 10* K, the threshold for 
atomic cooling. 

tions, and the resulting implications, are the subject of this 
Section. 

4.1 Computing properties of the dark matter 
hosts of the dwarf spheroidals 

To compute more rigorous estimates of properties of the 
dSphs' dark matter halos, we assume that the subhalo pop- 
ulation across all six Aquarius halos forms a representa- 
tive sample from ACDM simulations. We can then compute 
the distribution function of X (where X is, for example, 
Vmax, Vinfaii, Or Minfaii) for & dwarf by assigning a weight 
(likelihood jC) to each subhalo in our sample based on how 
closely it matches the measured M1/2 value of that dwarf. 
The posterior distribution of quantity X, given the measured 
value of Aii/2 and its error cta/, is given by 

P(XjMi/2, aw) oc P{X) £(X|Afi/2, aw) . (3) 

In practice, we compute moments in the distribution of 
quantity X via 

J2 X^P{X^)CiX,\A'h/2,^M) 

iX'^) = (4) 

J2 P{X,)C{X,\M^;2,'tm) 
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Table 2. Derived masses of the bright Milky Way dSphs. 
Columns: (1) Name of dwarf; (2) V-band luminosity (from 
Wolf et al. 2010); (3) median Vmax; (4) median Vinfaii; (5) me- 
dian Minfaii- The quantities in columns 3-5 are computed using 
Eqn. (4); the errors represent the 68.3% confidence interval. 



Name 


Lv 






Vinfall 


A'/infall 




[Lq] 




[kms-i] 


[kms"-"-] 






Fornax 


1 7-1-0.5 
-^•'-0.4 


X 


107 


i7.8«:^ 


22.0+4-1 


7.4+11 X 


iqs 


Leol 




X 


10® 






5.6l||? X 


iqs 


Sculpt 


9 c:+0.9 
^•^-0.7 


X 


10® 




iy + 5.8 

' -4.6 


6.61^:8 X 


iqs 


LeoII 




X 


10^ 


12.8+11 


le.oit; 


2.4+11 ^ 


iqs 


Sextans 




X 


10^ 








iqs 


Carina 




X 


10^ 


nA+_i:i 


14 4+^-"^ 


^■°-0.9 ^ 


iqs 


UrsMin 




X 


10^ 


2o.otl:t 


25 5+7-'' 




109 


CVnl 


9 q-l-0.4 
^••^-0.3 


X 


10^ 




14.51^.^ 


1.91?:° X 


iqs 


Draco 


9 9-1-0.7 
^•^-0.6 


X 


10^ 


20.5+_ti 


25 


i.2ig:?x 


109 



to compute properties of the hosts of the MW dwarf 
spheroidals. We assume that the hkelihood functions are log- 
normal, and compute the relevant moments in log as well, 

e.g., (logVmax) and (yfogV^^^ = {(logV^ax)^ - (logVniax)^)- 

The resulting values, with la errors, for Knax, 'Knfaii, 
and Mintaii are listed in Table 2. This analysis confirms 
the preliminary comparison performed in Section 3.1. The 
dwarfs are all consistent with subhalos in the range of 
10 < Vmax ^ 25kms~^, and at 95% confidence, none of the 
nine dwarfs have halos of Knax = 30kms^^ or greater. The 
central values of Vinfaii range from ~ 14 — 26 km s~^, with all 
but Draco and Ursa Minor having Vinfaii < 35kms~^ (2cr). 
Each simulation has 12-22 (16-33) subhalos with T^nfaii > 
35 (30) km s~^ and Knax > lOkms"^, however (the same is 
true for the Via Lactea II simulation: we find 27 subhalos 
with both Knax > lOkms"^ and Knfaii > 30kms^^). Even 
including for the Magellanic Clouds and Sagittarius, there 
are at least 6-21 subhalos in each simulation that are unac- 
counted for but have high enough masses that they should be 
luminous. Furthermore, the satellites predicted to be hosted 
by the most massive subhalos, Draco and Ursa Minor, are 
two of the three least luminous satellites in our sample. 

Eq. (4) assumes that the probability of a subhalo host- 
ing a specific dwarf depends only on how well that subhalo's 
(''1/2, dwarf) agrees with the measured M1/2 of the dwarf. 
If subhalos spanning a range of Knax or Knfaii have iden- 
tical values of M(ri/2, dwarf ), Eq. (4) assigns each of them 
the same probability of hosting that dwarf. In this case, the 
resulting mean value of Kmax may get larger weight from 
the numerous low-mass subhalos than from the less abun- 
dant high-mass subhalos since the mass function of CDM 
subhalos scales as A'^(> K) oc V~^. We have repeated our 
analysis with an additional weighting factor of to mimic 
a prior of equal probability per unit log V (effectively giving 
massive subhalos higher weights), and find that our results 
are qualitatively unchanged. This is not surprising, as the 
data strongly constrain M1/2, which, for subhalos, is tightly 
correlated with Kmax: the results are driven by the data, not 
by the choice of prior. 

Strigari et al. (2007b), Madau et al. (2008a), and 
Kuhlen (2010) have previously used the Via Lactea sim- 
ulations to derive constraints on Kmax (and, in the case of 



Kuhlen 2010, Knfaii) values. These were based on the masses 
of the satellites within either 300 or 600 pc (M300 or A/eoo)- 
Our calculations are based on a larger sample of ACDM 
subhalos and use more recently determined dynamical con- 
straints on the dwarfs - the masses with their de-projected 
half light radii - that have smaller errors than do M300 or 
Meoo- Furthermore, we have attempted to correct for the 
numerical effect of gravitational softening that affects M300 
at the ~ 20% level in simulations with the force resolu- 
tion of Aquarius level 2. This results in a decrease in the 
derived Knax values relative to the uncorrected case: ha- 
los of a fixed Knax have larger Ad^oo values after applying 
the correction. While our approach is somewhat more de- 
tailed than that of Madau et al. (2008a) and Kuhlen (2010), 
our results are generally consistent with those studies, and 
with Strigari et al. (2007b). Our results are also consistent 
with those of Strigari, Frenk, & White (2010), who have 
tested whether photometric and kinematic data on five of 
the bright dSphs are in accord with the gravitational po- 
tentials of Aquarius subhalos. They found good matches in 
each of the five cases, but always in systems with Kmax val- 
ues of 10 — 30kms~^, never in more massive halos. We find 
somewhat smaller Knax values for Draco, Leo I, Ursa Minor, 
and Sculptor than Peiiarrubia et al. (2008), likely because 
they adopted the Knax — ^max relation for field dark matter 
halos: Springel et al. (2008) find that the Aquarius subhalos 
are systematically denser than field halos, with an offset of 
0.2 dex in rmax at fixed Kmax- 

4.2 Comparison to ACDM predictions 

We are now in a position to directly compare observations 
and theoretical predictions for the hosts of the MW dSphs. 
The left-hand panel of Fig. 6 shows the luminosity func- 
tions of the Milky Way and M31 (solid black and gray 
curves) compared to the predicted luminosity functions from 
the Aquarius simulations using abundance matching, with 
Mi,/Lv ~ 2 (dotted curves). The abundance matching re- 
lation itself is computed by equating subhalo abundances 
from the Millennium and Millennium-II Simulations (closely 
following Guo et al. 2010) to galaxy abundances from the 
SDSS (Li & White 2009), with a power-law extrapolation 
to lower halo masses. While few of the simulated halos have 
subhalos massive enough to host the MW's brightest satel- 
lites, the Magellanic Clouds (Boylan-Kolchin et al. 2010, 
2011a; Busha et al. 2011), the agreement on the dSph scale 
{Mv > —14) is remarkably good in 5 of the 7 halos. Note 
that this agreement is not built into the abundance match- 
ing model: it results from the MW satellite luminosity func- 
tion having both the same slope and amplitude as the ex- 
trapolation of the field galaxy luminosity function to (much) 
lower luminosities than can currently be probed statistically, 
which is non-trivial. 

The masses of the subhalos that abundance matching 
predicts should host the dSphs are very different from those 
of the observed dSphs, however; this important point is illus- 
trated in the right-hand panel of Fig. 6. The black squares 
with errors show the Kmax values of the MW dSphs derived in 
Sec. 4.1, while the magenta circles show the measured Knax 
values from one of the numerical simulations that matches 
the luminosity function well (also colored magenta in the left 
panel of Fig. 6). The Kmax values of the simulated subhalos 
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Figure 6. Left: Observed luminosity functions for the Milky Way and M31 (thick solid lines) compared to abundance matching predictions 
based on the Aquarius simulations (thin lines, with Aq-E plotted in magenta; M-^jLy = 2 is assumed). Right: Values of Vmax computed 
in Sec. 4.1 for the nine luminous Milky Way dwarf spheroidals (square symbols with errors), along with Vmax(2 = 0) values of the 
subhalos with My < —8 (magnitudes are assigned by abundance matching) from the halo that best reproduces the luminosity function 
in the left panel (Aq-E). While numerical simulations combined with abundance matching reproduces the luminosity function of MW 
satellites, the structure of the dwarf spheroidals hosts' in this model does not match observations; the simulated subhalos are much more 
massive (have larger values of Vmax) than the dSphs. 



are systematically higher than those of the MW dSphs. It 
is therefore not possible to simultaneously match the abun- 
dance and structure of the MW dSphs in standard galaxy 
formation models based on dissipationless ACDM simula- 
tions. While there are many subhalos that match the struc- 
ture of the bright MW dSphs, these are not the subhalos 
that are predicted to host such galaxies in ACDM. 

The observed densities of MW satellites are very diffi- 
cult to reconcile with ACDM-based galaxy formation mod- 
els, where the stellar content of a galaxy is strongly cou- 
pled to l^infaii. To highlight the problem, we plot the in- 
ferred star formation efficiency ~ e^, = Af^/(/i, Minfaii), where 
fb = rib/^m is the universal baryon fraction - as a function 
of Minfaii in Fig. 7. The ellipses show 1 a uncertainties (note 
that the direction of the ellipses is due to the inverse cor- 
relation between and Minfaii at fixed M*). This relation 
is well-constrained at « = in the context of abundance 
matching for > 10*'^ Mq (approximately the complete- 
ness limit of the Li & White (2009) stellar mass function, 
corresponding to Mhaio = 6 x 10^" M©). The relation for 
Af* lower than the SDSS completeness limit is extrapolated 
using a power law (dashed portion of abundance matching 
lines) . 

The M* — Mhaio relation cannot be tested statistically 
on mass scales relevant for the dSphs at present, but it is 
immediately apparent that galaxy formation must proceed 
differently at Mhaio J; 10^" Mq than for larger systems if 
simulated subhalos accurately reflect the densities of the 
halos of dSphs as they exist the Universe. For example, the 
most luminous dSph of the MW, Fornax, has an inferred star 



formation efficiency of e* « 0.2, a value that is approached 
only at the scale of MW-mass halos. Ursa Minor and Draco, 
which are ~ 40 — 80 times less luminous than Fornax, sit 
in halos that are comparable or slightly more massive, and 
therefore have inferred efficiencies of closer to — 0.002. 



5 DISCUSSION 

Sections 3 and 4 have demonstrated that the structure and 
abundance of bright Milky Way satellites are not consis- 
tent with populating the most massive subhalos in hosts of 
Mvii « (1 — 2) X 10^^ Mq. In this Section, we discuss some 
possible remedies for this problem, ranging from downward 
revisions of the MW's dark matter halo mass (Sec. 5.1) to 
changes to ACDM (Sec. 5.4). 

5.1 Mass of the Milky Way 

The simulated halos used in this paper range from Mvir = 
9.5 X 10" to A/vir = 2.2 X 10^^ Mq. The true mass of the 
Milky Way is still a matter of significant uncertainty, how- 
ever. The apparent lack of massive subhalos might be under- 
standable if the Milky Way is significantly less massive than 
this simulated range. Ifere we summarize recent estimates 
of the Milky Way halo mass. 

• halo tracers 

Xuc et al. (2008) used blue horizontal-branch stars from 
the Sloan Digital Sky Survey, combined with mock obser- 
vations of hydrodynamical simulations of Milky Way-like 
galaxies, to find Mvir.Mw = l-Olo:2 x 10^^ M©, and M(< 
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Figure 7. The relation between Minf^n and star formation effi- 
ciency = Mi./(/(, Minfaii), the measure of the fraction of avail- 
able baryons converted into stars in a halo. Shaded ellipses show 
the nine bright dSphs considered in this paper, encompassing 
the 1 a uncertainties in M^^f^n from Table 2. The black shaded 
band and its extrapolation to lower masses (dashed lines) show 
e* (Afinfaii) derived from abundance matching. All of the bright 
dSphs have conversion efficiencies that are much higher than ex- 
pected from abundance matching, given their masses. Although 
the dSphs have similar values of Mi^faii, their inferred conver- 
sion efficiencies vary by two orders of magnitude. At larger mass 
scales, this spread is closer to a factor of two at fixed Minf^ii. 



60kpc) = (4.0 ± 0.7) X 10" Mq. Through a Jeans analysis 
of halo stars obtained from a survey for hypervelocity stars, 
Gnedin et al. (2010) found Mvir.MW = (1.6±0.3) x 10^^ Mq, 
and M{< 80kpc) = e.g!?;" x lO^M©. The largest uncer- 
tainties in these studies are the velocity anisotropy /3 and 
density profile (slope) assumed for the halo stars. Both Xue 
et al. (2008) and Gnedin et al. (2010) find most Ukely val- 
ues for P that are near 0.4, i.e., biased toward radial orbits. 
A smaller value of /3 would result in increases in the mass 
estimates. 

• satellite kinematics 

Kinematics of satellite galaxies in the Milky Way provide 
constraints on its mass. These constraints are also sensitive 
to the assumed velocity anisotropy of the satellite population 
and to whether or not Leo I, which has a large line-of-sight 
velocity, is considered a bona fide satellite. Watkins et al. 
(2010) find a "best estimate" - including Leo I - of M(< 
300 kpc) = (2.7 ±0.5) x 10^^ M©, which assumes that the 
satellite orbits are tangentially biased (based on observed 
proper motions); this estimate becomes (1.4±0.3) x 10^^ Mq 
for isotropic orbits, and (1.2 ± 0.2) x 10^^ Mq for radially 
anisotropic orbits (consistent with the Xue et al. and Gnedin 
et al. results, which also assume such orbits). 

• dynamics of the Large Magellanic Cloud 
Measurements of the proper motion of the LMC place it at 
a velocity of ~ 360kms~^ (Kallivayalil et al. 2006; Piatek 
et al. 2008), with statistical errors of less than 5%. This 



result implies that the Magellanic Clouds are likely on their 
first passage about the MW (Besla et al. 2007) , although this 
conclusion may depend on the mass of the MW (Shattow 
& Loeb 2009). Analyses of large cosmological simulations 
show that objects like the LMC are extremely unlikely to 
be unbound to their hosts (Boylan-Kolchin et al. 2011a), 
which allows for estimates of the MW's virial mass based 
on the dynamics of the Magellanic Clouds. Results of such 
estimates range from (1.2 — 2.0) x 10^^ Mq, depending on the 
treatment of baryonic physics and assumptions about the 
masses of the Clouds (Boylan-Kolchin et al. 2011a; Busha 
et al. 2011). 

• masers 

Masers can be used to constrain the angular speed of the 
local standard of rest, which, when combined with an esti- 
mate of the distance to the Galactic center, give the circular 
velocity of the MW at the sun's location. Reid et al. (2009) 
used this technique to derive V = 254 ± 16kms~^, signif- 
icantly higher than the lAU standard value of 220 km s"'^ 
(and therefore indicating a massive dark matter halo) . Ando 
et al. (2011) also used the maser method but found a circu- 
lar velocity of 213 ± 5kms~^, in strong disagreement with 
the Reid et al. value. This discrepancy is partially, but not 
completely, explained by different adopted distances to the 
Galactic Center _Ro: using i?o = 8.3 kpc (as in Reid et al.), 
Ando et al.'s estimate moves to 227kms~^. A combined 
Bayesian analysis by Bovy et al. (2009), which incorporates 
maser distances but uses greater freedom in the modeling 
than was allowed for in either Reid et al. or Ando et al., 
gives Vc(-Ro) = 236 ± llkms"^. This is higher than, but 
marginally consistent with, the lAU value. 

Converting between Vc(i?o) and Mvir is non-trivial. We 
can take some guidance from the MW mass models of 
Klypin, Zhao, & Somerville (2002), however. They found 
that the lowest halo mass that reasonably fit the MW was 
Mvir = 7.1 X 10" M0, which gave K(iio) = 216kms"^ 
(if angular momentum exchange between baryons and dark 
matter is included) or 246kms~^ (without angular momen- 
tum exchange) . These values require a maximal disk of = 
6 X 10^° Mq, and a value of M(< 100 kpc) = 3.8 x lO" Mq 
that is much lower than the more recent values of Gnedin 
et al. (2010) and McMillan (2011). The favored models of 
Klypin et al. (2002) assume Mvir = 1.0 x 10^^ M© and 
l/nrax « 220- 230kms-\ 

• other measures 

McMillan (2011) adopted a Bayesian approach to constrain- 
ing the mass of MW, incorporating photometric and kine- 
matic data as well as theoretical expectations from ACDM 
and modeling of the Galaxy. The result was a best-fitting 
mass within 100 kpc of M(< 100 kpc) = (0.84 ± 0.09) x 
lO^^M©, andavirial mass of Mvir = (1.26±0.24) x 10^^ M©. 
Based on abundance matching, Guo et al. (2010) found 
Mmw = 1.99 X 10^^ Mq, with an 80% confidence interval 
of (0.8 - 4.74) X 10^^ Mq. (These numbers are based on an 
overdensity criterion of A = 200; using our value of A = 94 
moves this range to (0.95 - 5.7) x 10^^ Mq.) 

The consensus of these results, then, is that the MW has 
a virial mass of ~ (1.0-2.0) X 10^^ Mq. The Aquarius halos 
used in this paper are broadly consistent with this range. 
Three of the halos - A, C, and D - lie at or slightly above 
2 X 10^^ Mq , while one - B - Ues just below 10^^ Mq . That all 
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of these halos have a substantial population of subhalos that 
are unaccounted for by the bright dSphs is a strong sign that 
the existence of "massive failures" is not an artifact of the 
choice of halo masses for the Aquarius halos. Nevertheless, it 
is important to consider the possibility that the MW's halo 
is less massive than is currently indicated by data. 

The substructure abundance in within the virial radii of 
dark matter halos is essentially self-similar when measured 
in terms of Vmax/K-ir (e.g., Gao et al. 2004; De Lucia et al. 
2004; van den Bosch et al. 2005; Gao et al. 2011), and the 
cumulative abundance of dark matter subhalos in a given 
halo scales as A''(> Knax) oc Vj^ax (Springel et al. 2008), so 
the number of halos of fixed Knax will scale linearly with 
Mvir. If the Aquarius halos are more massive than the halo 
of the Milky Way, then our first-order expectation is that the 
number of subhalos that are unaccounted for by the bright 
MW satellites will decrease by an amount proportional to 
the necessary reduction in virial mass. 

Since we find at least 6 subhalos with Vinfaii > 30 km s^^ 
that are unaccounted for when Mvir ~ 9.5 x 10^^ PIq and 16- 
21 such subhalos when Mvir = 2.2 x IO^^Mq, a reasonable 
scaling seems to be A^extra ~ 6 (Afvir/lO^^ M©). Reducing 
the halo mass to 5 x 10^^ Mq would then result in Aextra ~ 3, 
which perhaps could be explained by halo-to-halo fiuctua- 
tions. However, for the Milky Way to have a virial mass of 
5 X 10^^ Mq, all of the following would need to be true: (1) 
the LMC, SMC, and Leo I are all unbound to the Milky Way; 
(2) the Magellanic Clouds are extreme outliers in terms of 
their (high) masses; (3) stars in the Galactic halo are have 
high (radially biased) velocity anisotropics; (4) either the 
Local Group has a mass that is substantially smaller than 
what is derived based on the timing argument (Li & White 
2008), or M31 is 3-4 times more massive than the Milky Way 
in spite of having a fairly similar rotation velocity. 

Even if this is the case, it would still not explain the 
gap in Vmax between the most massive dwarf spheroidals (< 
25kms^^; Fig. 6, right-hand panel) and Magellanic Clouds 
(Vciic ^ 50kms~^): ACDM simulations do not contain such 
a gap. We therefore believe that a downward revision of 
the MW's mass is not the most likely explanation of the 
massive failures. As noted in Sec. 2.1, the abundance and 
structure of subhalos in the VL-II simulation (the only other 
simulation that has published data resolving subhalos on 
scales relevant to structure of the MW dwarfs) is very similar 
to that of the Aquarius simulations, providing no evidence 
that the background cosmology of the simulations causes the 
discrepancy. 

Finally, in the unlikely scenario where the Milky Way's 
mass has been systematically over-estimated, the unex- 
pected span in implied star formation efficiencies demon- 
strated in Fig. 7 would still persist. The only way to address 
the massive failures problem (Figs. 2-6) and the efficiency 
problem (Fig. 7) simultaneously is to modify the inner den- 
sities of massive dark matter subhalos significantly, and in 
a way that is anti-correlated with stellar mass. In the next 
section we explore whether baryonic feedback can plausibly 
remedy the situation. 



5.2 Baryonic feedback? 

Energy or momentum-driven feedback from star formation 
is expected to drive large-scale outfiows in galaxies at a va- 



riety of masses and redshifts. Although this feedback cannot 
couple hydrodynamically to dark matter at the centers of the 
halos of the Milky Way dwarfs, the removal of substantial 
amounts of gas would cause a dynamical re-arrangement, 
thereby reducing the central density of dark matter. The 
problems discussed above may be remedied if the central 
regions of the rotation curves of massive subhalos can be 
sufficiently reduced by such blow-out. Here, we aim to esti- 
mate the magnitude of this effect due to kinetic outflows. 

A supernova will release ~ 10^^ ergs of energy into its 
surroundings, a fraction esN of which will couple to the gas. 
If we parametrize the number of supernova explosions per 
100 solar masses of stars formed as A''ioo, then the gas mass 
carried by a supernova-driven outflow (with velocity Vout) 
due to episode of star formation is 



Afblow- 



4 A^ioo esN 



Vout 



500 km s" 



M* 



(5) 



The quantity in brackets is the mass-loading factor rj, and 
is expected to be of order = 1 — 2 for energy-driven winds 
and 77 oc Vrnax (though not exceeding ~ 5) for momentum- 
driven winds (e.g., Martin 1999; Springel & Hernquist 2003; 
Murray et al. 2005; Dave et al. 2011). 

A outflow occurring on a time scale that is short rela- 
tive to a halo's dynamical time effectively imparts an impul- 
sive change to energy of the dark matter particles on scales 
small relative to the characteristic size of the pre-outflow 
gas distribution. We model this effect by numerically com- 
puting an isotropic, spherically symmetric equilibrium dis- 
tribution function f{E) of a dark matter halo, modeled as a 
Hernquist (1990) profile, in the presence of a central gaseous 
component that is modeled as a generalized Hernquist pro- 
file (Dehnen 1993; Tremaine et al. 1994; see Boylan-Kolchin 
& Ma 2007 for a demonstration of the numerical stability 
of these models). The effect of impulsive blow-out of gas is 
modeled by removing the gaseous component and allowing 
the dark matter component to dynamically evolve to a new 
equilibrium in an A''-body simulation. This may be viewed 
as an instantaneous blow-out scenario, similar to that envi- 
sioned by Navarro et al. (1996) and Read & Gilmore (2005). 
As our experiments do not include any contraction of the 
dark matter halo in response to the initial assembly of the 
gaseous component, they should provide upper limits on the 
re-arrangement of dark matter in the blow-out model. We 
adopt a host halo of Vmax = 35kms~^ with rmax = 2kpc 
and investigate removing either 2.2 x 10^ or 1.1 x 10* Mq 
of gas. The scale radius of the gaseous component is set to 
200 pc, and the dark matter halos are resolved with 2 x 10® 
particles. 

Figure 8 shows the results of these experiments. Blow- 
out causes an immediate drop in the gravitational poten- 
tial at the system's center. Particles with r < i?s,gas lose a 
non-negligible amount of binding energy and move to larger 
apocenter distances on a dynamical time-scale, while par- 
ticles with pericenters Vp ^ /?s,gas are unaffected. The net 
effect is a reduction of the dark matter density for r < /?s,gas. 
Impulsively removing a gas mass of ~ 2 x 10^ Mq from 
the center of a typical Vinfaii ~ 35kms~'^ halo has a neg- 
ligible effect on the dark matter density profile (dashed 
line), while removing ~ 10* Mq of gas reduces the central 
Vcirc by approximately 25%. (The inner density cusp slope, 
which is originally —1, ends up closer to —0.5 after this 
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Figure 8. Simulations of impulsive blow-out from a Vmax = 
35kms~^ halo. The initial halo profile is plotted as a solid line; 
the dashed and dotted curves correspond to the final, relaxed 
profile after gas blow-out of 2.2 X 10^ and 1.1 X 10* Mq, respec- 
tively. Measured values of V^irc for the bright dSphs are plotted as 
squares, sized proportional to log Ly, with error bars. We empha- 
size that the two data points closest to the halo line post-blowout 
are among the least luminous dwarfs we consider (Draco and Ursa 
Minor, with Ly ~ 2 — 4 X 10^). Matching their densities via im- 
pulsive feedback would then require ejecting ~ 100 times as much 
mass as is present in stars today in these systems. 



maximal blow-out.) This result is consistent with Navarro 
et al. (1996), who used numerical simulations to show that 
a supernova-driven outflow following a single episode of star 
formation that created ~ 10* Mq of stars can create a core 
in the dark matter halo of a Vdrc ~ 50kms~^ dwarf galaxy. 
It also agrees with the work of Read & Gilmore (2005), 
which found that removing the vast majority of a dwarf 
galaxy's baryons (95-99%) in either one or two impulsive 
episodes significantly flattens a dark matter density cusp, 
with two blow-outs creating a dark matter core. Gnedin & 
Zhao (2002) demonstrated that removal of a gaseous disk 
containing nearly all of a dwarf galaxy's baryons will reduce 
the central dark matter density by approximately 50%; our 
reduction of ~ 25% in Vdic is in good agreement with this 
result. 

The final dark matter profile after removing 10* Mq of 
gas agrees only with Draco among the dSphs (square sym- 
bols with error bars). Blowing out 10* Mq of gas requires 
forming rj'^ 10* M© of stars, however [see Eqn. (5)]. Since 
the bright dSphs have M* « 5 x 10^ - 5 x 10^ Mq and r; 
should be of order unity (and certainly not in excess of 10), 
only the brightest dSph (Fornax) has formed enough stars to 
remove the necessary amount of gas. In particular, the nec- 
essary amount of star formation exceeds the stellar content 
of Draco by a factor of ~ 200, and the reduction in density 
is not sufficient to match Fornax's observed Vdic (the largest 
square, with Vdic ~ 17kms^^ at ~ 0.9kpc). 

Starting with a much more concentrated baryonic com- 
ponent - e.g., with a scale radius of 1-10 pc rather than 100 



pc, as assumed in this work - could cause gas blow-out to 
have a more deleterious effect on the central region of the 
dark matter halo. Collecting 10* Mq of gas in such a small 
region would likely steepen the initial dark matter profile a 
great deal, however. It would also result in star formation 
that is very centrally concentrated, in contrast to the typical 
half-light radii of the bright dwarfs of 300 < ^ 1000 pc. 
If the stars are indeed formed at very small radii, they would 
have to migrate outward in order to match observed sizes 
of the dSphs, which would perhaps result in radially biased 
velocities in the dwarfs' central regions. We conclude that 
absent outflows that are much larger than expected in usual 
models of galaxy formation (i.e., absent mass loading fac- 
tors of ?; = 10 — 100) or an extremely concentrated gaseous 
distribution, impulsive gas removal due to supernova feed- 
back is likely insufficient to lower the central densities of 
the massive subhalos predicted by ACDM enough to agree 
with observations of the bright dSphs. In the future, it may 
be possible to more directly constrain the amount of mass 
removed from dSphs through detailed analysis of chemical 
abundance patterns (Kirby et al. 2011). 

Governato et al. (2010) have recently performed a very 
high resolution simulation with a star formation density 
threshold of 100 cm~'^ and found it produces a cored galaxy 
with M^, = 5 X 10* Mq in a halo of Mvir = 5.7 x 10^° Mq 
{VciTc ~ 60kms~^). Pontzen & Governato (2011) argued 
that the effects of repeated cycles of star formation, blow- 
out, and re-cooling can explain these results. Unless this 
process is somehow much more efficient than the instan- 
taneous blow-out scenario discussed above - i.e., unless it 
removes much more gas mass per stellar mass formed than 
does a single episode of instantaneous blow-out - it is un- 
likely to explain the low densities of the Milky Way dwarf 
spheroidals: the Governato et al. galaxy has a stellar mass 
that is a factor of 10 larger than any of the MW dSphs (and a 
factor of 1000 larger than several of the bright dSphs). High 
resolution hydrodynamical simulations of the Milky Way's 
satellites find that the net effect of star formation and feed- 
back is either negligible or causes an increase in the dark 
matter density (di Cintio et al. 2011; Parry et al. 2012). 
Clearly, it will be invaluable to have a larger number of sim- 
ulations of possible dSphs over a wider mass range, simu- 
lated with different codes and stellar feedback prescriptions. 
It will also be interesting to better constrain the possible 
role of black hole feedback, which can power strong outfiows 
without relying on star formation, in the Milky Way dSphs. 
Reines et al. (2011) have recently discovered an AGN in a 
dwarf starbursting galaxy with a mass similar to that of the 
LMC (i.e., more massive than the MW dSphs); as of yet, 
there is no evidence for central black holes in the dSphs, 
with Jardel & Gebhardt (2012) finding a 1 a upper limit of 
3.2 X 10" Mq for the mass of a potential black hole at the 
center of Fornax. 



5.3 Stochastic galaxy formation? 

If baryonic feedback has not strongly modified the struc- 
ture of massive MW subhalos, and the abundance of these 
objects is commensurate with that found in the Aquar- 
ius simulations, then it seems unavoidable that galaxy for- 
mation must be highly stochastic in halos of < 50kms~'^ 
(Afinfaii < 10^" Mq). By this, we mean that the stellar mass 
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of a galaxy cannot correlate with its dark matter halo mass 
at these scales (see Fig. 7). Stochasticity is not unexpected 
in low halo masses, as gas cooling will depend strongly on 
the ability to form molecular gas, which in turn depends 
on the gas metallicity (Kuhlen et al. 2011). An additional 
source of scatter in stellar mass at fixed halo mass comes 
from allowing satellites to have a Mhaio — Af* relation that 
differs from that of central galaxies (Neistein et al. 2011). 
Yang et al. (2011) have demonstrated that such a variation 
is expected because the stellar mass of a satellite should be 
related to the Mhaio — M* relation at Zinfaii, not at z = 0. (It 
remains to be seen whether the abundance of MW satellites 
is correctly reproduced in this or similar abundance match- 
ing models.) 

However, the stochasticity required by the results of 
Table 2 and Fig. 7 is somewhat more curious: it requires a 
strong (and systematic) suppression of galaxy formation in 
the most dense subhalos of the Milky Way, which is counter- 
intuitive. It also requires that the luminosities of all subhalos 
with 30 < Vmax ^ 60 km s~^ - a range in which we expect to 
find several objects in the Milky Way - host galaxies with 
luminosities a factor of 10-1000 lower than subhalos in the 
range 15 < Vmax ^ 30kms~^. If this is due to scatter, it 
must be only scatter in one direction, toward low M* at 
fixed Afhaio, in the mass range of massive failures, which 
would host galaxies with luminosities comparable to ultra- 
faint dwarf spheroidal galaxies in this scenario. 

The dynamics of the ultra-faint are not as well con- 
strained as those of the bright dSphs for two main rea- 
sons. First, the kinematic data for the ultra- faint galax- 
ies of the Milky Way are not of the same quality as for 
the bright dSphs. Additionally, the ultra-faints have smaller 
half-light radii, which means that measurements of Mi/ 2 
have less power in constraining Knax. It is therefore possi- 
ble that some, or all, of the massive failures are accounted 
for by ultra-faint galaxies; at least, this is not presently 
ruled out by the data. If this is the case, however, then 
galaxy formation must be extremely stochastic in halos of 
Vinfaii J; 40 — 50kms~^, as it would require M/L varia- 
tions in excess of 1000 at this scale to simultaneously ex- 
plain the ultra-faints and the brightest dSphs. As shown in 
Fig. 4, the massive failures are all more massive than the 
photo-suppression threshold at all times, so it is not clear 
what shuts off star formation in the massive failures hosting 
ultra-faints in this scenario. For example, the hybrid simu- 
lations of Bovill & Ricotti (2011a,b), which reproduce many 
observed properties of the ultra-faint dwarf spheroidals by 
placing them in halos below the filtering mass, over-produce 
brighter Milky Way dwarfs, even though only star formation 
prior to reionization is considered. 

5.4 Different dark matter physics on sub-galactic 
scales? 

Should none of the explanations in Sections 5.1 - 5.3 prove 
correct, then modifications of the ACDM model on small 
scales may be necessary to explain the structure of the MW 
dSphs. Indeed, Lovell et al. (2012) have already examined 
the structure of a version of the Aquarius A halo simulated 
in a Warm Dark Matter cosmology. They find that WDM 
substantially alleviates the massive failure problem, as the 
subhalos in this simulation are (1) less numerous and (2) less 



dense than their counterparts in the CDM run. Although 
this result is very promising, it remains to be seen whether 
WDM simulations produce enough substructure to account 
for the ultra-faint galaxies of the MW (Maccio & Fontanot 
2010; Polisensky & Ricotti 2011). The particle mass used by 
Lovell et al. (2012) is also near the lower limit of what is al- 
lowed by observations of Lyman-a absorption features in the 
spectra of distant quasars (Boyarsky et al. 2009). Improved 
measurements of this type will further constrain the allowed 
parameter space for WDM and the viability of AWDM as a 
cosmological model. 

Alternatively, it may be the case that dark matter is 
not strictly collisionless, but rather has a non-negligible 
cross section for self-scattering (Spergel & Steinhardt 2000). 
This self-interacting dark matter (SIDM) model was orig- 
inally proposed to resolve the missing satellites problem 
and to explain the apparently overdense hosts of low sur- 
face brightness galaxies predicted by CDM, but interest in 
SIDM mostly abated when it was shown that SIDM would 
produce strong effects in clusters of galaxies (among other 
reasons; see, e.g., Yoshida et al. 2000). Nevertheless, the 
range of models explored was relatively restricted, and Loeb 
& Weiner (2011) have recently argued that SIDM with a 
velocity-dependent cross section peaking at a velocity scale 
relevant for dwarfs retains the benefits of original SIDM 
models while avoids their drawbacks (see also Feng et al. 
2009). Simulating structure formation in such a model would 
therefore be of significant interest. 



5.5 Missing Physics at V < 50kms ^? 

Ultimately, the structure of the Milky Way dwarfs must 
be understood as part of a full theory of galaxy formation 
that can explain the formation, evolution, and properties of 
all galaxies. It is therefore reasonable to ask whether the 
structural issues we have raised regarding MW dSphs are 
also seen in other systems of similar mass. Although ob- 
serving galaxies at these masses over cosmological volumes 
is essentially impossible at present, THINGS (Walter et al. 
2008) has observed 34 nearby galaxies {D < 15 Mpc) in HI 
with high spatial resolution. We use the results of Oh et al. 
(2011a,b), who performed a detailed mass modeling of seven 
low-mass THINGS galaxies, to compare the Mi, — Vinfaii re- 
lation for these isolated galaxies to the bright MW dSphs 
and to abundance matching expectations in Fig. 9. 

Intriguingly, the THINGS data mostly lie above the 
abundance matching relation and its extrapolation to lower 
masses and provide a smooth transition from the abun- 
dance matching M* — Vinfaii relation at Vinfaii ^ 75kms~'^ 
to the regime of the MW dSphs at ~ 15 — 25 km s~^ derived 
in this work. Simulations of isolated dwarf galaxies with 
masses comparable to the THINGS sample also fall above 
the abundance matching curve (Sawala et al. 2011). Note, 
however, that if the ACDM model is correct, and galajc- 
ies with l^nfaii ~ 40kms~^ typically have stellar masses of 
~ 10® Mq , then current stellar mass functions are substan- 
tially underestimating the abundance of such galaxies. Al- 
ternatively, either most halos with Vinfaii ~ 40kms~^ host 
systems with much lower total stellar masses, or the abun- 
dance of halos at these masses is lower than what is predicted 
by ACDM. There are at least three additional lines of evi- 
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Figure 9. Inferred Mt(V'infaii) relation for bright MW dSphs 
(black squares with error bars, based on calculations in Sec. 4.1). 
The Magellanic Clouds (right-pointing triangles) are placed on 
the plot at their current values of Vflat, which is a lower limit 
to Vinfall- Observations of low-mass field galaxies from THINGS 
(tabulated in Oh et al. 2011a) are plotted as open black stars. 
These galaxies all lie higher than the z = abundance match- 
ing relation (solid curve), as well as its extrapolation to lower 
Mnfall (dashed curve), and the deviations are systematically larger 
at lower values of Vinfall- The shaded region around the abun- 
dance matching relation shows a scatter of 0.2 dex in M* at 
fixed Vinfaii, which is the upper limit allowed for massive halos 
Mnfall > 150kms~i (Guo ct al. 2010; Bchroozi ct al. 2010). 



dence arguing that isolated halos with Vd 
not match ACDM expectations: 



50kms do 



• HI observations: The ALFALFA survey has per- 
formed a blind 21-cm emission line search over a wide area 
to look for neutral hydrogen in galaxies. Papastergis et al. 
(2011) have shown that while the velocity width function 
measured from ALFALFA agrees fairly well with ACDM pre- 
dictions for massive galaxies, the observed number counts 
fall below those predicted by ACDM for w < lOOkms"^ 
(corresponding approximately to Vmax ^ 75kms~^, assum- 
ing an average conversion of w^q — 0.75 Vmax)- The discrep- 
ancy reaches a factor of ~ 8 at w = 50kms~^ (Vmax ~ 
37kms~^) and becomes even worse at lower w. 

• Void galaxies: Tikhonov & Klypin (2009) analyzed 
properties of voids in the Local Volume in comparison to 
theoretical predictions. They find that the abundance of void 
galaxies is over-predicted by a factor of ~ 10 in ACDM at 
Vcirc ~ 40kms~^, and that the void size distribution is only 
reproduced if halos of Vcirc ^ 40kms~^ do not host void 
galaxies (but see Tinker & Conroy 2009). 

• Damped Lyman-a systems: The gaseous content of 
dark matter halos at 2 ~ 3 can be probed by quasar absorp- 
tion spectra. Barnes & Haehnelt (2009) have shown that 
many of the properties of damped Lyman-a systems can be 
understood in ACDM-based models. This success comes at 
the expense of requiring halos with T4irc 50kms~^ to be 



very baryon-poor. As noted in Section 3.3, 50kms~^ is well 
above the photo-suppression scale at this redshift, indicat- 
ing that reionization should not have caused such halos to 
lose a substantial amount of their baryons. 

There are, of course, many potential sources of these dis- 
agreements, and the underlying ACDM theory is certainly 
not the most likely of of these sources. A better understand- 
ing of feedback from star formation and its effects on halos 
of Vinfaii 50kms~^ will be crucial, and may explain all 
of these apparent discrepancies, as well as other issues such 
as the central densities of low surface brightness galaxies 
(though see Kuzio de Naray & Spekkens 2011 for arguments 
against baryonic physics explaining the density structure of 
these galaxies). It is imperative not to rely on plausibility 
arguments for the effects of feedback, but rather to under- 
stand whether realistic feedback models can actually pro- 
duce dwarf spheroidal galaxies with properties akin to those 
seen in the Milky Way (as challenging as this may be!). 



6 CONCLUSIONS 

In this paper, we have expanded on the arguments of BBK, 
where we first showed that the bright satellites of the Milky 
Way apparently inhabit dark matter subhalos that are sub- 
stantially less dense that the most massive subhalos from 
state-of-the-art ACDM simulations. Using subhalo profiles 
computed directly from the simulations rather than assum- 
ing subhalos are fit by NFW profiles, we have confirmed 
our previous result. Furthermore, we have now computed 
the most likely Knax, Vinfaii, and Minfaii values of the dwarf 
spheroidals using a likelihood analysis of the Aquarius data. 
This procedure predicts that all of the Milky Way dwarf 
spheroidals reside in halos with Knax J; 25kms~'^, whereas 
more than ten subhalos per host halo are expected to have 
Knax > 25kms~^. 

This "massive failure" problem cannot be solved by 
placing the bright satellites in the subhalos with the largest 
values of Knax at infall or at the epoch of reionization: as dis- 
cussed in Section 3.3 and shown in Figures 4-5, the missing 
subhalos are among the most massive at all previous epochs 
as well. Explaining this lack of galaxies in the expected mas- 
sive subhalos is not natural in standard ACDM-based galaxy 
formation models: options include (1) a Milky Way halo that 
either is significantly deficient in massive subhalos, or is pop- 
ulated by subhalos with much lower concentrations than are 
typical; (2) stochasticity in galaxy formation at low masses, 
such that halo mass and luminosity have essentially no corre- 
lation; (3) strong baryonic feedback that reduces the central 
density of all massive subhalo by a large amount (> 50% 
reduction on scales of ~ 0.5 kpc). 

We have argued above that these solutions all seem 
fairly unlikely as individual causes. It might be possible to 
apply them all at once: if the Milky Way halo mass sits at 
the low end of current constraints (~ lO'^^ Mq), and galaxy 
formation produces order unity scatter in M* at fixed halo 
mass below ~ 50kms~^, and baryonic feedback is able to 
alter the central densities of dark matter halos in a maxi- 
mal way, then it may be possible to explain the low den- 
sities of the MW dSph galaxies. We find this combination 
somewhat implausible, but it is certainly worth exploring. 
A detailed comparison of the masses of M31 dSph galaxies 
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will be particularly useful in assessing the " rare Milky Way" 
hypothesis. Finally, if we do reject solutions (l)-(3), then we 
are left with the question of the nature of dark matter. Al- 
lowing for phenomenology such as self-interactions, decays, 
or non-negligible thermal velocities may explain the puzzles 
discussed here without destroying the many successes of cold 
dark matter models on large scales. 
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APPENDIX A: MASS PROFILES OF 
SIMULATED SUBHALOS 

In BBK, we used the values of Knax and rmax computed 
in Springel et al. (2008) to construct NFW profiles for the 
Aquarius subhalos. While subhalos are generally well-fitted 
by NFW profiles in the aggregate, any individual subhalo 
may show (significant) deviations from the NFW profile 
specified by its {Vmax, f'max} values. In this paper, there- 
fore, we use the raw particle data to compute the circular 
velocity profiles of the Aquarius subhalos. 

While using the raw particle data provides a much bet- 
ter estimate of the true mass profile of a subhalo than does 
the NFW assumption, we still must remember that we are 
interested in properties of subhalos on scales that are not 
much larger than the force resolution of the Aquarius sim- 
ulations. Draco, for example, has a de-projected half-light 
radius of 291 pc (Wolf et al. 2010), while the Plummer- 
equivalent softening length e of the Aquarius simulations is 
65 pc. The effect of force softening is to reduce the density 
on scales of ~ 3 e or smaller (note that forces are Newtonian 
on scales larger than 2.8 e). This will not affect differential 
quantities such as p(r) on larger scales, but will affect cu- 
mulative properties such as Vdvdr), such that the measured 
mass is an underestimate of the true mass. For example. 
Font et al. (2011) find that the mass within 300 pc of the 
centers of Aquarius A-2 subhalos (the resolution level of the 
simulations used in this paper) is systematically underesti- 
mated by approximately 20% by comparing with the higher 
resolution Aquarius A-1 simulation. 

We adopt the following procedure to correct for the ef- 
fects of force softening. First, we fit the density profile of each 
subhalo with an Einasto profile on the range^ [291 pc, ruppor] , 
where rupper is the smaller of 3 kpc or 1.5 rmax- This proce- 
dure typically results in very good fits: using the goodness- 



^ We choose 291 pc for the lower limit of our fits because this is 
the de-projected half-light radius of Draco, the most dense dSph 
in our bright sample. 



Figure Al. The difference between the measured subhalo density 
profiles and the best-fitting Einasto density profiles as a function 
of radius. The mean values (symbols), along with the standard 
deviations (error bars) based on the ten highest Vmax subhalos 
in the Aquarius B halo are plotted. The subhalos follow Einasto 
profiles (all with goodness-of-fit Q values of < 0.03), with very 
small residual scatter, for r > 250 pc. On smaller scales, however, 
the measured profile falls below the Einasto fit because of grav- 
itational softening. This causes an underestimate of ycirc('') for 
r < 2.8 € (forces for r > 2.8 e are exactly Newtonian); this value 
is marked with a vertical dashed line. 

of- fit measure 

1 '^"'"^ 

= -^r i^OgPsnbsiu) - log/9modcl(ri)]^ , (Al) 

we find a median (mean) Q value of 0.057 (see also Fig. Al). 
We then use the best-fitting Einasto profile to model the 
density distribution for r < 291 pc and the raw particle dis- 
tribution for r > 291 pc. The correction to Vdic from this 
procedure will be of decreasing importance with increasing 
r; the correction is typically 5-10% at 291 pc. For a small 
number of subhalos (typically with low numbers of parti- 
cles), the best-fitting Einasto profile has Q > 0.1. In this 
case, we opt to be conservative and do not apply any cor- 
rection to the measured particle distribution. We also do not 
apply the correction for subhalos with rmax < 500 pc. 

To assess the impact of this softening correction on our 
results, we re-ran the analysis from Section 4.1 without in- 
cluding the correction; Table Al compares the two sets of 
results. The softening correction has minimal impact on the 
derived values of Knax, with none of the dSphs changing by 
more than 5%. Vinfaii and Minfaii are very stable as well: 
only Draco, the dSph with both the smallest half-light ra- 
dius and largest inferred Vmax, is affected at the 9% level in 
Knfaii (20% in Minfaii). Eveu these changes are much smaller 
than the errors on these values. We therefore conclude that 
our results are not driven by including (or excluding) the 
correction for gravitational softening. 
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Table Al. Effects of the softening correction on derived prop- 
erties of the Milky Way dwarfs. For each dSph, we list two sets 
of numbers: the first row includes the softening correction, while 
the second row does not include the softening correction. 
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vmax 


''mfall 


Mini 


Fall 




[kms-1] 


[kms-1] 


[Mq] 




Fornax 
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17.9l«f, 
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Leo I 


— 2a) 
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